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We  describe  a  new  algorithm  for  computing  the  efficient  frontier  of  the  "bi-objective  maximum-flow 
network-interdiction  problem."  In  this  problem,  an  "interdictor"  seeks  to  interdict  (destroy)  a  set  of  arcs  in 
a  capacitated  network  that  are  Pareto-optimal  with  respect  to  two  objectives,  minimizing  total  interdiction  cost 
and  minimizing  maximum  flow.  The  algorithm  identifies  these  solutions  through  a  sequence  of  single-objective 
problems  solved  using  Lagrangian  relaxation  and  a  specialized  branch-and-bound  algorithm.  The  Lagrangian 
problems  are  simply  max-flow  min-cut  problems,  while  the  branch-and-bound  procedure  partially  enumer¬ 
ates  s-t  cuts.  Computational  tests  reveal  the  new  algorithm  to  be  one  to  two  orders  of  magnitude  faster  than 
an  algorithm  that  replaces  the  specialized  branch-and-bound  algorithm  with  a  standard  integer-programming 
solver. 
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1.  Introduction 

In  the  deterministic  maximum-flow  network-inter- 
diction  problem  (MXFI),  an  "interdictor"  wishes  to 
restrict  an  "adversary's"  use  of  a  capacitated  net¬ 
work.  Specifically,  the  adversary  seeks  to  maximize 
flow  through  the  network,  while  the  interdictor 
seeks  to  minimize  that  maximum  flow  by  inter¬ 
dicting  (destroying)  arcs  using  available  interdiction 
resources  (Wood  1993).  These  resources  could  be  a  set 
of  cruise  missiles,  a  given  number  of  aerial  sorties, 
etc.  MXFI  and  related  problems  appear  in  the  areas 
of  military  planning  (Whiteman  1999),  drug  interdic¬ 
tion  (Steinrauf  1991),  and  protecting  civil  infrastruc¬ 
ture  against  terrorist  attacks  (Salmeron  et  al.  2004). 
The  model  may  be  classified  as  a  bilevel,  target- 
selection,  or  weapons-allocation  model;  see  Matlin 
(1970),  Bracken  and  McGill  (1974),  Bracken  et  al. 
(1977),  and  Wen  and  Hsu  (1991).  A  simpler  version 
of  MXFI,  the  problem  of  cutting  all  flow  using  mini¬ 
mum  interdiction  resource,  provided  one  of  the  earli¬ 
est  applications  of  maximum  flows  and  the  max-flow 
min-cut  theorem.  The  flow  was  rail  traffic  from  the 
Soviet  Union  into  eastern  Europe;  see  Harris  and  Ross 
(1955)  and  Ford  and  Fulkerson  (1956,  1957). 

This  paper  studies  an  extension  of  MXFI,  which 
we  call  the  "bi-objective  maximum-flow  network- 
interdiction  problem"  (BMXFI).  In  BMXFI,  the  inter¬ 
dictor  wishes  to  identify  all  Pareto-optimal  solutions 
for  MXFI,  i.e.,  the  efficient  frontier,  with  respect  to 
minimizing  post-interdiction  maximum  flow  and 


minimizing  "total  interdiction  cost,"  i.e.,  interdiction- 
resource  expenditure.  A  bi-objective  formulation  is 
important  for  a  military  commander  who  (i)  must 
pre-plan  for  various  resource  availabilities  given  the 
uncertainty  of  warfare,  or  (ii)  may  wish  to  consider 
cost-to-effectiveness  or  risk-to-effectiveness  tradeoffs 
that  are  difficult  to  include  directly  in  an  interdiction¬ 
planning  model. 

Steinrauf  (1991)  and  Wood  (1993)  show  that  MXFI 
can  be  formulated  as  an  integer  program  (IP).  Thus, 
BMXFI  can  be  solved  as  a  set  of  IPs  for  MXFI,  one 
for  each  plausible  level  of  interdiction  resource.  But, 
MXFI  is  strongly  NP-complete  and  difficult  to  solve 
in  practice  (Wood  1993).  Hence,  solution  times  could 
be  unacceptable,  especially  in  a  time-constrained,  mil¬ 
itary  setting.  Benders  decomposition  can  be  helpful, 
because  it  incorporates  easy-to-solve  maximum-flow 
subproblems,  but  the  computational  expense  becomes 
prohibitive  because  so  many  subproblems  must  be 
solved  (Cormican  1995).  Derbes  (1997),  Bingol  (2001) 
and  Uygun  (2002)  devise  fast  heuristics  based  on 
Lagrangian  relaxation,  but  these  cannot  guarantee 
good  solutions.  Implicitly,  these  authors  also  consider 
BMXFI,  and  construct  heuristic  solution  approaches 
based  on  weighted-sums  scalarization  of  the  two 
objective  functions. 

The  solution  of  a  sequence  of  IPs  will  typically  re¬ 
quire  that  a  commercial  solver  be  installed  on  one 
or  more  local  computers,  and  that  up-to-date  licenses 
be  maintained  for  those  computers.  But,  analysts  at 
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a  hastily  organized  military  headquarters  should  not 
need  to  worry  about  software  licenses.  Alternatively, 
an  open-source  IP  solver  could  be  used,  but  such 
solvers  tend  to  be  less  efficient  than  their  commercial 
counterparts.  Consequently,  a  need  exists  for  a  spe¬ 
cialized  approach  to  solving  BMXFI. 

This  paper  constructs  a  procedure  for  solving 
BMXFI  based  on  weighted-sums  scalarization  of  the 
objectives,  viewed  as  an  application  of  Lagrangian 
relaxation.  Lagrangian  relaxation  partially  maps  out 
the  efficient  frontier  through  a  sequence  of  max- 
flow  min-cut  problems.  The  procedure  then  fills 
in  the  missing  pieces  of  the  frontier  with  a  spe¬ 
cialized  branch-and-bound  algorithm  that  enumer¬ 
ates  near-minimum-capacity  cuts  with  respect  to 
Lagrangianized  capacities.  For  simplicity,  we  phrase 
our  discussion  in  terms  of  "Pareto-optimal  solutions" 
and  "the  efficient  frontier."  However,  for  a  fixed  limit 
on  total  interdiction  cost,  we  typically  allow  a  1%  or 
5%  optimality  gap.  Thus,  we  are  actually  identifying 
"near-Pareto-optimal  solutions"  (Carlyle  et  al.  2003), 

1. e.,  the  "nearly  efficient  frontier." 

The  next  section  defines  MXFI  and  BMXFI  precisely. 
Sections  3  and  4  show  how  Lagrangian  relaxation  and 
cut  enumeration  can  solve  MXFI  for  specific  limits 
on  total  interdiction  cost.  Section  5  presents  the  full 
solution  procedure.  Section  6  presents  computational 
examples,  and  Section  7  presents  conclusions. 

2.  Problem  Formulation 

Let  G  —  (N,  A)  denote  a  directed  graph  with  node  set 
N  and  arc  set  A  c  N  x  N.  Two  special  nodes  s  f  t  are 
identified:  s  is  the  source  node  and  t  is  the  sink  node. 
We  also  define  FS(i)  —  {(/',  f)  e  A\i'  —  i)  and  RS(i )  = 
[{j',i')zA\i'  =  i}. 

Each  arc  keA  has  a  capacity  uk  e  2+  and  interdic¬ 
tion  cost  rk  e  %+.  The  interdiction  cost  is  the  amount  of 
some  resource  necessary  to  destroy  arc  k,  i.e.,  reduce 
that  arc's  capacity  to  zero.  We  assume  that  rk  >  0  and 
uk  >  0  for  all  ke  A,  since  any  arc  with  rk  =  0  or  uk  —  0 
may  be  viewed  as  "already  interdicted"  and  can  be 
removed  from  the  problem.  We  also  expect  that,  in  a 
military  setting,  all  rk  will  be  small  integers.  For  exam¬ 
ple,  they  might  represent  the  number  of  aerial  sorties, 
one,  two  or  three,  say,  required  to  ensure  destruction 
of  a  bridge.  We  note  that  only  simple  variations  are 
required  to  model  an  undirected  network,  partial  arc 
interdictions,  and  node  interdiction  (Wood  1993). 

An  s-t  cut  is  a  partition  (Ns,Nt)  of  N  such  that 
s  eNs  and  t  eN,.  All  cuts  in  this  paper  are  s-t  cuts, 
so  we  drop  "s-t"  hereafter.  Given  a  cut  C  =  (Ns,  Nt), 
Ac  denotes  the  set  of  arcs  k  —  ( i ,  j)  e  A  such  that  i  e  Ns 
and  j  eNt.  The  capacity  of  C,  defined  by  J2keAc  uk>  Pro~ 
vides  an  upper  bound  on  the  maximum  s-t  flow  in  G; 
deletion  of  Ac  from  G  disconnects  all  paths  directed 


from  s  to  t  (e.g.,  Ahuja  et  al.  1993,  p.  185).  A  cut  C  is 
minimal  if  Ac  is  a  minimal  disconnecting  set.  We  use 
%  to  denote  the  collection  of  all  minimal  cuts  in  G. 

Let  xk  —  1  if  arc  k  is  interdicted,  and  let  xk  —  0, 
otherwise.  Also,  let  G(x)  =  (N,  A  —  {k  e  A  \  xk  =  1}), 
and  let  g(x)  denote  the  maximum  s-t  flow  in  G(x).  In 
BMXFI,  an  interdiction  plan  x  e  {0, 1}|A|  is  dominated  if 
there  exists  another  plan  x'  such  that  (i)  y(x')  <  g(x), 
(ii)  J2keA  rkx'k  —  X keA  rkxk>  and  (id)  at  least  one  inequal¬ 
ity  is  strict.  An  interdiction  plan  that  is  not  dominated 
is  Pareto-optimal.  We  define  BMXFI  to  be  the  problem 
of  identifying  all  Pareto-optimal  interdiction  plans. 

Let  yk  denote  the  flow  on  arc  k  e  A,  and  let  y„ 
denote  the  flow  on  an  artificial  "return  arc"  a  —  (f,  s) 
f  A.  Also,  define  the  standard  shorthand  notation 
"minxeX{y(x),  /i(x) }"  to  mean:  Find  all  Pareto-optimal 
pairs  (g(x),  h(x))r  for  xeX,  with  respect  to  minimiza¬ 
tion.  We  can  then  formulate  BMXFI  based  on  the  well- 
known  maximum-flow  problem: 


BMXFI 


min 

xe{0, 1)IAI 


(1) 


where 


g(x)=  max  ya  (2) 

y  s£»W,y, 

s-t.  8,ya  +  £  yk-  £  y*  =  0  VieN  (3) 

fceFS(f)  JfceRS(i) 

0<xjk<uk(l-xk)  VkeA,  (4) 

with  8S  —  —1,  8t  —  1,  and  S,  =  0  for  all  i  e  N  —  {s,  t). 
Equations  (2)-(4)  define  a  standard  maximum-flow 
problem  with  "interdiction-modified  arc  capacities," 
uk(  1  —  xk),  for  all  k  e  A. 

We  identify  Pareto-optimal  interdiction  plans  by 
isolating  one  objective  function  and  limiting  the  value 
of  the  other  by  a  constraint.  Isolating  the  max¬ 
imum-flow  objective  and  imposing  a  limit  R  on 
total  interdiction  cost,  we  obtain  the  single-objective 
maximum-flow  network-interdiction  problem: 


MXFI(R) 

z*(R)  =  ming(x),  (5) 

xeX 

where  X=|xe  {0, 1}|A|  |  rkxk  <  R}. 

MXFI(R)  can  be  converted  to  a  simple,  minimizing 
IP  by  reformulating  the  problem  into  this  essentially 
equivalent  model  (Cormican  et  al.  1998) 


z*  (R)  =  min  max  ya  -  £  xkyk 

X£X  y  keA 

s.t.  (3)  and  0  <yk<uk, 


(6) 

VkeA,  (7) 
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and  by  then  taking  the  dual  of  the  inner  problem. 
(See  Wood  1993,  also.)  An  optimal  0-1  solution  always 
exists,  so  the  resulting  IP,  after  a  few  simplifications,  is 

MXFI-IP(R) 

z*(R)  =  min  J2ukPk 

s.t.  al-aj  +  l3k  +  xk>0  V(i,j)  =  keA  (8) 

E  rk*k  <  R  (9) 

keA 

all  variables  e  {0, 1) 

as  =  0,  at  =  1. 

The  full  efficient  frontier  can  be  identified  by  solv¬ 
ing  MXFI-IP(R),  using  standard  IP  software,  over  a 
sufficiently  wide  range  of  "R-values."  However,  this 
approach  can  be  computationally  costly,  so  we  present 
an  alternative  approach  in  Section  5.  This  approach  is 
also  based  on  solving  MXFI(R),  so  we  first  focus  on 
this  single-objective  problem. 

Another  useful  view  of  MXFI(R)  derives  from  the 
max-flow  min-cut  theorem  (Ford  and  Fulkerson  1956): 
Replace  the  max-flow  problem  that  defines  g(x)  with 
an  equivalent  min-cut  problem.  This  leads  to  the 
following  proposition,  in  which,  without  loss  of  gen¬ 
erality,  we  have  interchanged  the  order  of  the  two 
minimizations. 

Proposition  1.  MXFI(R)  can  be  solved  by  finding 
C*  e  %  such  that  a  maximum-capacity,  cost-feasible  set  of 
arc  interdictions  in  C*  leaves  as  little  uninterdicted  capac¬ 
ity  as  possible.  That  is,  MXFI(R)  is  equivalent  to 

MXFI-C(R) 

z*(R)  —  min  min  E  uk0  ~xk)  (10) 

(_£©  X 

ksAc 

s.t.  E  rkxk  <  R  (11) 

keAc 

xk  e  {0, 1}  Vke  Ac  (12) 

xk  —  0  Vk  e  A-  Ac.  (13) 

Proof.  As  in  the  dual  to  the  maximum-flow  prob¬ 
lem,  any  feasible  solution  to  MXFI-IP(R)  identifies  a 
minimal  or  nonminimal  cut  C  =  (Ns,Nt)  by  setting 
a,  —  0  for  all  i  e  Ns,  and  a,  =  1  for  all  i  e  N,  (e.g..  Wood 
1993).  Constraints  (8)  are  satisfied  by  setting  (3k  —  1  or 
xk  =  1  (but  not  both)  for  all  k  e  Ac,  and  by  setting  (3k  = 
xk  —  0  for  all  ke  A  —  Ac.  Thus,  if  we  extend  So  in  (10) 
to  include  nonminimal  cuts,  MXFI-C(R)  and  MXFI- 
IP(R),  and  thus  MXFI(R),  are  equivalent.  The  validity 
of  restricting  So  to  minimal  cuts  is  obvious.  □ 

For  a  given  cut  C,  the  inner  minimization  of  MXFI- 
C(R)  comprises  a  knapsack  problem  with,  for  our  pur¬ 
poses,  small  coefficients  rk  e  %+.  Thus,  if  the  number 


of  minimal  cuts  were  small  (unlikely  in  practice), 
MXFI(R)  would  be  fairly  easy  to  solve:  Enumerate 
all  minimal  cuts  (e.g..  Shier  and  Whited  1986),  and 
solve  one  knapsack  problem  for  each.  If  rk  =  1  for  all 
k  e  A,  the  problem  may  appear  to  be  easier,  because 
the  knapsack  problem  becomes  trivial,  but  it  remains 
strongly  NP-complete  (Wood  1993). 

3.  Lagrangian  Relaxation 

We  approach  the  solution  of  MXFI(R)  by  moving 
the  interdiction-cost  constraint  (11)  into  the  objective 
function  using  Lagrangian  relaxation  (Derbes  1997). 
Let  A  denote  the  Lagrangian  multiplier  associated 
with  (11).  Then,  for  any  R  e  j£+  and  A  e  [0,  oo),  we 
define 

MXFI-LR(A,  R) 

z(A,  R)  =  minmin  Eu/c(l  —  xk) 

keAc 

+  a(Ew-k)  (14) 

\keAc  / 

s.t.  xke{0,l}  Vke  Ac  (15) 

xk  =  0  Wk  e  A  -  Ac.  (16) 

This  Lagrangian  problem  yields  a  lower  bound  for 
MXFI(R):  For  any  R  e  %+  and  A  e  [0,  oo),  z(A,  R)  < 
z*(R).  Note  that  if  a  solution  (C,  x)  to  MXFI-LR(A,  R) 
satisfies  Y.keA*  rA  —  Rr  then  z(A,  R)  =  z*(R),  that  is, 
x  is  optimal  for  MXFI(R). 

Let  /(A)  denote  the  max-flow  value  in  G  given 
arc  capacities  u'k  —  m  i  n  j  uk ,  \rk\.  Observing  that  (i)  the 
objective  function  in  MXFI-LR(A,  R)  may  also  be 
written  as  J2k€Ac(uk  +  (Ah  —  uk)xk)  —  AR,  and  (ii)  the 
problem's  sole  constraints  are  (15)  and  (16),  a  simple 
solution  appears  for  the  inner  minimization  of  MXFI- 
LR(A,  R):  For  each  k  e  Ac,  set  xk  —  Hi  A rk  <  up,  other¬ 
wise  set  xk  =  0.  Hence, 

z(A,  R)  =  min(  minjfij.,  \rk)  —  AR  )  (17) 

Ce<®  \kSAc  / 

—  /(A)  —  AR,  (18) 

where  the  last  equality  follows  from  the  max-flow 
min-cut  theorem.  Thus,  evaluating  z(A,  R)  is  no 
harder  than  solving  a  max-flow  problem.  We  do  not 
use  the  following  fact,  but  it  is  interesting  to  note 
that  the  best  Lagrangian  bound,  maxA>0  z(A,  R),  can 
be  computed  in  polynomial  time;  see  Lawler  (1976), 
pp.  94-97,  Megiddo  (1979),  and  Derbes  (1997). 

For  brevity  in  the  following  sections,  let  the  inner 
minimizations  in  MXFI-C(R)  and  MXFI-LR(A,  R)  be 
defined,  respectively,  as 

z*(R,  C)  =  min  J2uk0-~xk)  (19) 

X  keAc 

s-t.  E  Vk  <  R 

keAc 
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xke{0,l}  VkeA 
xk  =  0  Wk  e  A  —  Ac, 

and 

z(A,  R,  C)  =  minjMfc,  Xrk}  —  AR.  (20) 

keAc 

4.  Cut  Enumeration 

We  will  not  need  to  maximize  the  Lagrangian  lower 
bound  z(A,  R)  for  every  R  when  solving  BMXFI. 
However,  it  is  useful  for  now  to  imagine  that,  for  a 
given  R,  we  have  identified  A  >  0  such  that  z(A,  R) 
is  "reasonably  close"  to  maxA>0  z(A,  R).  In  the  process 
of  identifying  A,  we  will  have  identified  a  cut  C  that 
establishes  that  bound,  i.e.,  z(A,  R)  =  z(A,  R,  C),  and 
will  have  found  a  feasible  solution  x  to  MXFI(R)  that 
defines  an  upper  bound  z(R)  >  z*(R).  In  particular, 
having  identified  C,  we  can  set  z(R)  —  z*(R,  C). 

Now,  if  z(R)  —  z(A,  R)  <  e  for  some  prespecified  tol¬ 
erance  e  >  0,  we  have  identified  an  e-optimal  solution 
for  MXFT(R).  If  not,  we  can  apply  Theorem  1,  below, 
to  find  one.  We  first  establish  two  lemmas  to  simplify 
the  theorem's  proof. 

Lemma  1.  z(A,  R,  C)  <  z*(R ,  C)  for  any  A  e  [0,  oo), 
R  6  Z+,  and  C  e  %. 

Proof.  This  follows  because  (i)  z*(R,  C)  is  the  opti¬ 
mal  objective  value  for  MXFT(R)  defined  on  a  network 
that  consists  only  of  the  arcs  Ac  connected  in  parallel 
from  s  to  f,  and  (ii)  z(A,  R,  C)  defines  a  Lagrangian 
lower  bound  on  the  same  problem.  □ 

Lemma  2.  Let  A  e  [0,  oo)  and  Re%+.  If  z(R)  is  a  con¬ 
stant  such  that  z(R)  >  z*(R),  then  the  outer  minimiza¬ 
tion  in  MXFI-C(R )  can  be  restricted  to  C  e  'S  such  that 
z(A,  R,  C)  <  z(R)  without  changing  that  problem's  opti¬ 
mal  value. 

Proof.  Suppose  that  Ce?  is  a  solution  of  the  outer 
minimization  in  MXFT-C(R)  with  z(A,  R,  C)  >  z(R). 
Then,  by  Lemma  1, 

z(R)  <  z(A,  R,  C)  <  z*(R,  C)  =  z*(R),  (21) 

which  is  a  contradiction.  □ 

Theorem  1.  Given  optimality  tolerance  e  >  0,  A  e 
[0,  oo),  and  R  e  %+ ,  suppose  that  Ce^  yields  upper 
bound  z(R)  —  z*(R,  C).  Also,  define 

z*(R)  =  min  z*(R,C),  (22) 

Ce«(z,A,R,s) 

where  ^(z.  A,  R,  e)  =  { C  e  %  \  z(A,  R,  C)  <  z(R)  —  e)  U 
{Cj.  Then,  z*(R)  —z*(R)  <  e,  i.e.,  z*(R)  is  an  e-optimal 
objective  value  for  MXFI(R). 

Proof.  If  z*(R)  >  z(R)  —  e  —  z*(R,  C)  —  e,  then  C 
has  already  yielded  an  e-optimal  solution  and  the 


theorem  is  valid.  Suppose  not,  i.e.,  suppose  z*(R)  < 
z*(R,  C)  —  e.  Then,  z*(R,  C)  —  e  is  a  valid  upper  bound 
on  z*(R)  and  thus  Lemma  2,  with  z(R)  replaced  by 
z*(R,  C)  —  e,  implies  that  the  minimization  in  (22) 
yields  z*(R).  □ 

Corollary  1.  If  C  is  a  minimizer  of  (22),  then  every 
minimizer  xl  o/(19),  with  C  replaced  by  C,  is  an  e-optimal 
interdiction  plan  for  MXFI(R).  □ 

When  the  Lagrangian  procedure  of  the  previous 
section  leaves  z(R)  —  z(A,  R)  >  e,  Theorem  1  and 
Corollary  1  guide  us  to  a  solution  of  MXFI(R): 
Using  A  =  A  that  approximately  maximizes  z(A,  R), 

(i)  enumerate  all  cuts  C  having  z(A,  R,  C)  <  z(R)  —  e, 

(ii)  solve  a  knapsack  problem  for  each  C  in  order 
to  minimize  uninterdicted  capacity,  i.e.,  compute 
z*(R,  C),  and  (iii)  save  the  solution  with  the  least 
uninterdicted  capacity.  Naturally,  whenever  a  better 
solution  is  found,  the  upper  bound  can  be  updated, 
and,  in  effect  then,  we  may  apply  Theorem  1  multi¬ 
ple  times  while  searching  for  an  e-optimal  solution. 
(Thus,  beginning  with  a  solution  that  is  not  e-optimal, 
and  a  correspondingly  weak  upper  bound,  does  not 
imply  that  we  must  then  seek  an  optimal  solution, 
which  the  theorem  might  seem  to  suggest.) 

We  implicitly  enumerate  all  necessary  cuts  using 
the  tree-search  algorithm  in  Balcioglu  and  Wood 
(2003).  That  algorithm  has  this  key  feature:  Given  that 
each  arc  k  e  A  has  capacity  u'k  (u'k  =  minji^,  A rk},  in 
our  case),  then  each  node  in  the  search  tree  corre¬ 
sponds  to  a  unique  cut  C  such  that 

(23) 

k&Q 

for  every  "child"  cut  C  of  C.  Our  use  of  this  algorithm 
defines  a  depth-first  branch-and-bound  algorithm. 
In  contrast  to  LP-based  (linear-programming-based) 
branch  and  bound  that  reoptimizes  the  lower  bound 
for  each  child,  our  algorithm  recomputes  a  valid 
bound  for  each  child,  but  does  not  reoptimize  with 
respect  to  A:  We  have  not  found  that  reoptimiza¬ 
tion  improves  computation  times,  in  practice.  Enu¬ 
meration  strategies  other  than  depth  first  could  be 
employed,  but  this  strategy  is  easy  to  program  and 
works  well  in  practice. 

We  outline  the  cut-enumeration  algorithm  next, 
using  | It]  to  denote  the  smallest  integer  at  least  as 
large  as  the  real  number  h.  This  ceiling  function  ap¬ 
pears  because,  given  our  integral  data,  a  nonintegral 
lower  bound  can  always  be  rounded  up  to  the  nearest 
integer. 

Procedure  Enumerate.  /*  Solves  M XFT ( R )  for 
specified  Re%+.  *  / 

Input  Data.  Data  for  MXFT(R),  feasible  solution  x(R) 
and  corresponding  upper  bound  z  (R),  parameters 
A  >  0  and  e  >  0,  and  global  lower  bound  [z(A,  R)~\ . 
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Output.  An  e-optimal  interdiction  plan  x(R)  for 
MXFI(R). 

0.  Set  arc  capacities  u'k  •<—  minjiij.,  A rk)  for  all  k  e  A; 

1 .  Begin  the  tree-search,  cut-enumeration  algorithm 
(Balcioglu  and  Wood  2003); 

2.  For  each  cut  C  encountered, 

(a)  If  z(R)  -  [z(A,  R,  C)]  =  z(R)  -  [(Zk£Ac  K  - 
AR)"|  <  e,  go  to  Next; 

(b)  Compute  z*(R,  C)  and  its  optimizer  x*c; 

(c)  If  z*(R,  C)  <  z(R),  set  z(R)  •<-  z*(R,  C)  and 
x{R)^x*c; 

(d)  If  z(R)  —  |"z(A,  R,  C)]  <  e  go  to  Next; 

(e)  For  each  child  cut  C  of  C, 

(i)  Continue  the  tree  search  recursively  at 
Step  2  with  C  <—  C; 

(ii)  Upon  returning  from  exploring  C,  if  z(R)  — 
|~z(A,  R,  C)]  <  e  go  to  Next; 

/*  The  above  check  may  be  useful  if  z(R)  has 
improved.  */ 

(f)  Next:  If  C  is  the  initial  cut  in  the  tree,  or  if 
z(R)-Tz(A,R)l  <  e,  halt  the  enumeration  and  Return 
x(R)  and  z (R); 

(g)  Backtrack  from  C. 

We  note  that  Procedure  Enumerate  may  actually 
identify  some  nonminimal  cuts  that  must  be  searched 
recursively,  but  which  are  not  candidates  for  yield¬ 
ing  optimal  solutions.  Our  implementation  adjusts  for 
this,  but  we  ignore  the  issue  above  for  simplicity. 


5.  A  Complete  Algorithm 

We  now  return  to  BMXFI.  The  full  efficient  frontier 
for  BMXFI  can  be  explored  by  solving  MXFI(R)  for 
all  possible  R-values,  for  example  by  repeated  appli¬ 
cation  of  Procedure  Enumerate.  We  construct  a  more 
efficient  approach,  however,  which  exploits  the  rela¬ 
tionship  between  MXFI(R)  and  MXFI(R'),  for  R  close 
to  R'. 


5.1.  Overview 

Using  the  arguments  in  the  proof  of  Proposition  1,  we 
find  that  BMXFI  is  equivalent  to 

BMXFI-C 


min  min 

Ce%  x 


E  «*(!-**)/  E  rk*k 

keAc  keAc 


s.t.  xke{0,l}  VkeAc 
xk  —  0  Wk  e  A  —  Ac. 


(24) 

(25) 

(26) 


In  view  of  this  and  ignoring  the  constant  A R,  we 
see  that  MXFI-LR(A,  R)  is  the  reduction  of  BMXFI- 
C  by  means  of  weighted-sums  scalarization  of  its 
objectives,  using  (unsealed)  weights  1  and  A.  Hence, 
the  evaluation  of  /(A)  (see  (18))  for  a  range  of  A-values 


is  equivalent  to  solving  weighted-sums  scalarizations 
using  various  weights. 

While  we  cannot  expect  to  identify  the  full  efficient 
frontier  this  way  (Climaco  et  al.  1997),  we  hope  that 
a  substantial  portion  can  be.  Thus,  for  various  values 
of  A,  we  perform  the  following  steps: 

1.  Using  a  maximum-flow  algorithm,  compute 

/(A)  min  E  minjiq.,  A rk]  and 

C€<B  k€Ac 

C  arg min  E  min{ uk,  \rk}; 

CeFS  keAc 

2.  Set  xk  —  1  whenever  uk  >  A rk  and  k  e  Aq,  and  set 
xk  *—  0  otherwise; 

3.  Set  R  4-  Y.keA  rA,  z*(R)  <-  T.keA  *4(1  -  xk),  and 
x(R)  < —  x. 

The  solution  x(R)  is  optimal  for  MXFI(R)  since 
YlkcA  rkxk(R)  —  R-  For  problems  MXFI(R)  not  immedi¬ 
ately  solved  using  these  steps,  we  have  a  head  start 
because,  for  any  R,  z(A,  R)  =  /(A)  —  AR  <  z*(R)  for  all 
values  of  A  that  the  procedure  examines. 

Let  Rmax  be  the  smallest  R  such  that  z *(R)  =  0,  i.e., 
the  minimum  amount  of  interdiction  resource  that  can 
force  the  adversary's  maximum  flow  to  0.  The  cor¬ 
responding  interdiction  plan  is  found  by  setting  A  = 
minkeA  uk/rk  in  the  previously  described  procedure. 
Also,  for  any  A  >  maxkeA  uk/rk,  the  solution  obtained 
corresponds  to  R  =  0.  Thus,  it  suffices  to  examine  val¬ 
ues  of  A  in  the  range  [min^  uk/rkl  A  +  maxfeA  uk/rk ] 
for  some  A  >  0. 

Since  /(A)  is  a  piecewise-linear  function  with  no 
more  than  Rmax  break  points,  a  procedure  akin  to 
Benders  decomposition  could  map  out  z(A,  R)  com¬ 
pletely  by  solving  at  most  Rmax  +  1  max-flow  prob¬ 
lems.  The  Benders  master  problem  would  be  simple 
enough  to  solve  by  inspection.  However,  the  sequence 
of  max-flow  min-cut  problems  would  involve  widely 
varying  arc  capacities  and  would  be  harder  to  solve 
efficiently  than  one  in  which  arc  capacities  are  non¬ 
decreasing.  Consequently,  we  opt  for  mapping  out 
the  function  /(A)  approximately,  starting  with  a  small 
value  for  A  and  increasing  it  according  to  an  empiri¬ 
cally  derived  rule. 

We  motivate  that  rule  with  a  network  consisting  of 
a  set  of  parallel  arcs  of  the  form  k  —  (s,  f),  all  with 
unique  ratios  uk/rk.  In  this  case,  /(A)  breaks  only  at 
points  where  some  arc  capacity  minjdj.,  A rk\  switches 
from  A rk  to  uk,  i.e.,  at  A  =  uk/rk.  Therefore,  we  can  map 
out  the  function  /(A),  and  thus  the  functions  z(A,  R) 
for  all  R,  by  evaluating  the  maximum  flow,  or  mini¬ 
mum  cut  capacity,  at  only  these  values  of  A. 

Naturally,  the  situation  is  more  complex  for  a  gen¬ 
eral  network.  But,  for  simplicity,  let  us  continue  to 
assume  that  the  ratios  uk/rk  are  unique.  Suppose 
some  A'  >  0  has  been  defined,  arc  capacities  u'k  — 


180 


Royset  and  Wood:  Solving  the  Bi-Objective  Maximum-Blow  Network-Interdiction  Problem 
INFORMS  Journal  on  Computing  19(2),  pp.  175-184,  ©2007  INFORMS 


min{uk,  \'rk]  have  been  computed,  and  a  minimum 
cut  C  is  identified  such  that  u'k  —  uk  for  at  least  one 
k  e  Ac.  Without  loss  of  generality,  assume  A 'rk  ^  uk 
for  any  k  e  A.  Now,  our  example  implies  that  /(A) 
must  break  somewhere  on  the  interval  (A',  A"]  where 
A"  =  min ksAc{uk/rk  \  uk/rk  >  A'},  i.e.,  between  its  cur¬ 
rent  value  and  the  value  where  the  function  would 
break  if  C  were  the  only  cut  in  the  network.  Since 
G  typically  contains  a  huge  number  of  cuts,  some 
other  cut  may  become  "the  minimum  cut"  before  A 
reaches  A".  Thus,  increasing  A  all  the  way  to  A"  tends 
to  overshoot  the  next  breakpoint.  This  observation, 
along  with  empirical  testing,  leads  us  to  adopt  the 
following  rule  for  computing  the  next  A: 

A  =  max{(l  —  +  4>\" ,  A'  +  Aj  (27) 

for  some  empirically  determined  A  >  0  and  <l>  e  (0, 1]. 
We  find  that  A  =  0.1  mint:gA  min{«A.,  rk },  and  </>  =  0.25 
work  well  in  general. 

5.2.  Algorithm  Interdict 

Based  on  the  discussion  and  outline  above,  we  now 
state  a  complete  algorithm  for  BMXFI: 

Algorithm  Interdict.  /*  Solves  BMXFI.  */ 

Parameters.  Absolute  tolerance  e  >  0  and  algorithm 
controls  <fi  e  (0, 1]  and  A  >  0. 

Input  Data.  Network  G  =  (At,  A),  with  interdiction 
costs  rk  e  %+  and  arc  capacities  uk  e  2+  for  all 
ke  A.  Source  node  seN  and  sink  node  t  e  At. 
Output.  An  optimal  or  e-optimal  solution  x(R)  to 
MXFI(R)  for  each  R  e  (0, 1, . . . ,  Rmax). 

0.  Set  A  4—  min keA{uk/rk}  and  ©4-0; 

/*  ©  represents  the  set  of  R-values  with  known 
solutions.  */ 

1.  While  A  <  maxteA{ uk/rk]  +  A 

/*  Step  1  roughly  maps  out  /(A)  and  identifies,  we 
hope,  a  large  number  of  solutions  x(R)  for  different 
R-values.  *  / 

(a)  Solve  the  maximum  flow  problem  on  G  using 
arc  capacities  u'k  =  minju,.,  Ar(. },  and  find  a  corre¬ 
sponding  minimum-capacity  cut  C; 

(b)  Set  xk  4-  1  for  all  k  e  Ac  such  that  uk  >  A rk, 
and  set  xk  <-  0  otherwise; 

(c)  Set  R  4-  J2keA  W 

(d)  If  R  ^  ©,  /*  A  new  solution  has  been  found.  */ 
set  ©  4-  ©  U  {R},  A(R)  4-  A,  z*(R)  4-  J2ksA  uk{l  -  xk), 
and  x(R)  4-  x; 

(e)  If  A  =  min keA{uk/rk},  set  Rmax  4-  R; 

(f)  If  minji/j.,  Xrk)  =  uk  for  all  k  e  Ac,  break 
from  "While"  loop;  Else,  set  A  4-  max{(l  —  <l>)\  + 
(l>mmk€Ac{uk/rk  \  uk/rk  >  A),  A  +  A);  /*  See  (27)  */ 

2.  For  each  R  e  {0, 1, ... ,  Rmax)  —  © 

/*  Step  2  refines  lower  bounds  for  R-values  whose 
e-optimal  solutions  remain  unidentified.  It  may  also 
identify  some  new,  absolutely  optimal  solutions.  */ 


(a)  Select  R',  R"  e  ©,  R'  <  R  <  R",  to  bracket  R  as 
tightly  as  possible; 

(b)  Set  A (R)  4-  argmaxA6{A(R/)  A(R„)1  z( A,  R); 

(c)  Compute  a  new  A-value,  A'"  4-  (z*(R')  — 
z*(R"))/(R"  —  R'); 

/*  A'"  estimates  argmaxA>0  z(A,  R)  using  subgra¬ 
dients  of  z(A,  R).  * / 

(d)  If  A'"  ^  (A(R"),  A(R'))  continue  the  "For"  loop; 

/*  The  "continue"  means  that  A'"  cannot  yield  an 

improved  lower  bound,  so  just  proceed  to  the  next 
R-value.  */ 

(e)  Compute  C  and  x  by  performing  Steps  1(a)- 
(b),  with  A  =  A'",  and  set  z(A"',  R)  4-  z(A'",  R,  C); 

(f)  Set  A(R)  4-  argmaxA6(A(R)/A,„)  z(A,  R); 

(g)  If  R  =  J2keA  rkxk, 

/*  A  new  optimal  solution  has  been  found.  */ 
set  ©  4-  ©  U  {Rj,  z*(R)  4-  z(A(R),  R),  and  x(R)  -4-  x; 

3.  For  each  R  e  {0, 1, ... ,  Rmax)  —  ©  in  increasing 
order  /*  Step  3  identifies  e-optimal  solutions  for 
R-values  that  have  eluded  exact  solution  in  Steps  1 
and  2.  Depending  on  e,  it  may  or  may  not  be  neces¬ 
sary  to  actually  call  the  cut-enumeration  (branch-and- 
bound)  procedure.  */ 

(a)  Set  initial  values  z(R)  4-  z(R  —  1),  x(R)  4- 
x(R  —  1),  A  4—  A (R);  /*  z(A,  R)  is  now  also  defined  * / 

(b)  If  z(R)  —  [z(A,  R)]  <  e,  then  continue  the 
"For"  loop; 

(c)  Call  Procedure  Enumerate  with  inputs  z (R), 
x(R),  A,  [z(A,R)],  to  return  e-optimal  x(R); 

4.  /*  Output  solution.  */ 

(a)  For  all  Re©,  Print(R,  "Optimal  solution  is," 

x(R)); 

(b)  For  all  R  e  {0, 1, . . . ,  Rmax}  —  ©,  Print(R, 
"e-optimal  solution  is,"  x(R)). 

Step  2  of  Algorithm  Interdict  estimates  the  best 
lower  bound  maxA>0z(A,R)  for  z*(R)  and  the  corre¬ 
sponding  maximizer  A(R)  by  using  z*(R'),  z*(R")  and 
a  slope  estimate  of  z(A,  R)  with  respect  to  A,  where 
R'  and  R"  bracket  R.  In  particular.  Step  2(c)  uses 
the  fact,  from  (18),  that  z(A,  R)  <  z(A,  p)  +  (p  —  R)A 
for  any  p  e  2+  and  A  >  0,  as  well  as  the  fact  that 
z(A(p),  p)  =  z*(p)  for  p  =  R',  R".  Step  2  may  not  max¬ 
imize  the  lower  bound  precisely,  but  optimizing  the 
bound  has  not  proved  computationally  worthwhile. 
We  also  note  that  Steps  1  and  2  are  actually  inte¬ 
grated  in  the  algorithm's  implementation,  and  have 
been  separated  above  for  the  sake  of  clarity. 

Step  3(a)  in  Algorithm  Interdict  exploits  the  fact 
that  an  upper  bound  on  z*(R')  is  also  an  upper  bound 
on  z*(R)  for  R  >  R'.  Flence,  the  initial  upper  bound 
for  each  R  is  computed  without  significant  computa¬ 
tional  cost. 

Since  Rmax  is  finite  and  the  network  G  has  a  finite 
number  of  cuts.  Algorithm  Interdict  solves  BMXFI  in 
finite  computing  time.  Given  that  the  number  of  cuts 
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can  be  exponential  in  the  size  of  G,  the  algorithm's 
worst-case  complexity  is  exponential,  however. 

5.3.  Enhancements  to  the  Algorithm 

Our  implementation  of  Algorithm  Interdict  makes 
two  modifications  to  the  basic  algorithm  to  improve 
computational  efficiency.  First,  maximum  flows  are 
calculated  for  many  similar  problems,  so  we  exploit 
"warm  starts."  We  use  a  variant  of  the  shortest- 
augmenting-path  algorithm  of  Edmonds  and  Karp 
(1972)  to  solve  maximum-flow  problems,  and  an 
inherent  feature  of  such  algorithms  is  that  the  maxi¬ 
mum  flow  for  a  problem  with  arc  capacities  uk]  pro¬ 
vides  a  feasible  initial  flow,  i.e.,  a  warm  start,  for  a 
problem  with  capacities  42)  >  for  all  k  e  A.  We 
make  use  of  that  feature  in  Step  1  of  the  algorithm 
since  arc  capacities  never  decrease  as  A  increases.  We 
also  make  use  of  it  in  Step  2  of  Procedure  Enumerate, 
where  the  next  cut  to  consider  is  determined  through 
a  warm-started  max-flow  calculation  (Balcioglu  and 
Wood  2003). 

The  second  computational  improvement  comes 
from  perturbing  the  ratios  uk/rk  so  that  uk/rk  ^  ut/rl 
for  all  k,l  e  A,  k  ^  1.  Derbes  (1997)  first  notes  the 
benefits  of  this  technique;  see  also  Bingol  (2001).  We 
implement  perturbations  by  multiplying  every  uk  by 
a  large  positive  integer  (in  practice,  6  x  105),  and  by 
then  adding  a  small  positive  integer  to  uk  when¬ 
ever  uk/rk  =  U\tr\  for  some  l  e  A,  l  ^  k.  The  pertur¬ 
bations  are  several  orders  of  magnitude  smaller  than 
the  smallest  arc  capacity  and  hence  do  not  change 
a  problem's  solutions.  This  technique  typically  cre¬ 
ates  more  linear  segments  in  z(A,  R)  as  a  function 
of  A.  Hence,  it  is  more  likely  that  z(A,  R)  is  constant 
for  A  in  some  interval.  Because  any  A  in  this  interval 
yields  z(A,  R)  =  z*(R),  this  technique  results  in  more 
instances  of  MXFI(R)  being  solved  optimally  in  Step  2 
of  the  Algorithm. 

6.  Computational  Results 

This  section  examines  the  efficiency  of  Algorithm 
Interdict  for  solving  BMXFI  applied  to  artificial  and 
real-world  network  data.  The  algorithm  is  coded  in 
Java  1.4.2  (Sun  Microsystems  2004),  and  run  on  a 
3.0  GHz  Pentium  IV  laptop  computer  with  one  giga¬ 
byte  of  RAM. 

We  compare  our  algorithm's  run  times  to  the  times 
required  by  a  hybrid  algorithm  that  looks  just  like 
Algorithm  Interdict  except  that  a  commercial  IP  solver 
is  called  instead  of  Procedure  Enumerate  to  resolve 
instances  of  MXFI-IP(R)  for  "problematic"  values 
of  R.  We  use  CPLEX  Version  9.1  (ILOG  2005)  as  the 
solver,  integrated  with  our  Java  code  using  ILOG 
Concert  Technology  (ILOG  2006).  We  have  tested  a 
variety  of  CPLEX  options,  and  find  that  it  is  best 


to  favor  "branching  up"  and  to  apply  mixed-integer 
rounding  cuts  aggressively. 

The  tables  of  results  also  report  the  average  number 
of  branch-and-bound  nodes,  computed  over  all  R  e 
{0, 1, . . . ,  Rmax),  for  both  algorithms.  When  CPLEX 
solves  an  IP  without  any  enumeration,  it  records  zero 
nodes  enumerated.  For  consistency,  we  also  record 
zero  nodes  when  Algorithm  Interdict  does  not  invoke 
Procedure  Enumerate  for  a  given  R. 

The  CPLEX  version  of  MXFI-IP(R)  exploits  z(R  —  1) 
just  as  Procedure  Enumerate  does,  and  incorporates 
two  enhancements.  The  first  defines  an  auxiliary  inte¬ 
ger  variable  x  and  an  auxiliary  "branching  constraint" 
J2keAxk  —  x  (Appleget  and  Wood  2000).  When  not  all 
rk  =  1,  setting  the  branching  priority  highest  for  x, 
(and  second  highest  for  the  xk,  and  lowest  for  the 
other  variables)  can  substantially  reduce  enumeration. 
Roughly  speaking,  this  causes  the  solver  to  branch  on 
the  number  of  interdictions  before  branching  on  indi¬ 
vidual  interdictions.  The  second  enhancement  incor¬ 
porates  warm  starts  for  the  sequence  of  closely  related 
LP  relaxations  that  CPLEX  must  solve.  We  believe 
that  these  enhancements  to  the  IP-based  methodol¬ 
ogy  mean  that  the  solution  times  reported  for  the 
hybrid  algorithm  present  a  substantial  challenge  to 
Algorithm  Interdict. 

Figure  1  illustrates  an  instance  from  the  first  class 
of  test  networks,  rectangular  grid  networks.  Letting 
nx  and  n2  denote  the  number  of  rows  and  columns 
of  nodes  in  these  networks,  \N\  —  nxn2  +  2  and  |A|  = 
2 +  2 ((«2  —  1  )nx  +  (n,  —  1  )n2).  Arcs  connecting  the 
source  s  and  the  sink  t  all  have  effectively  infinite 
capacities  and  are  invulnerable  to  interdiction.  All 
other  arcs  k  have  capacities  uk  that  are  drawn  ran¬ 
domly  and  independently  from  the  discrete  uniform 
distribution  on  [1,49]. 

We  consider  three  variants  of  the  grid  networks,  all 
having  common  topology  but  with  different  numeri¬ 
cal  data.  For  variant  Al,  interdiction  costs  are  rk  =  1 
for  all  k  e  A.  Interdiction  costs  are  randomly  gener¬ 
ated  for  variants  A2  and  A3,  as  specified  in  Table  1. 
In  this  table,  kw,  kE,  ks,  and  kN  represent  west¬ 
bound,  eastbound,  southbound,  and  northbound  arcs, 
respectively.  Randomly  generated  data  are  statisti¬ 
cally  independent.  Small  integers  should  adequately 
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Table  1  Probabilities  for  Pseudo-Random  Generation 
of  Interdictions  Costs  in  Grid-Network  Variants 
A1,A2,  and  A3 


Network 

Al 

A2 

A3 

II 

-O 

o 

ol 

1 

1/4 

1/2 

Prob(r%=2) 

0 

3/4 

1/2 

ProbK  =  3) 

0 

0 

0 

I! 

U. 

-C2 

o 

ol 

1 

0 

0 

CM 

II 

-Q 

o 

ol 

0 

1 

1/2 

CO 

II 

-Q 

o 

ol 

0 

0 

1/2 

ii 

_Q 

O 

Q_ 

1 

1/4 

1/2 

Prob(4s  =  2) 

0 

3/4 

1/2 

Prob(4s  =  3) 

0 

0 

0 

ProWkN  =  1) 

1 

1/4 

1/2 

Prob(r%=2) 

0 

3/4 

1/2 

Prob  (r%  =  3) 

0 

0 

0 

model  interdiction  costs  in  military  applications,  so 
rk  e  {1, 2, 3}  in  all  problem  instances.  Note  that  Algo¬ 
rithm  Interdict  can  be  simplified  when  all  rk  =  1 
because  the  knapsack  problem  in  Step  2(b)  of  Proce¬ 
dure  Enumerate  becomes  trivial.  We  have  not  imple¬ 
mented  code  to  take  advantage  of  this  special  case, 
however. 

Tables  2,  3,  and  4  display  solution  times  and  aver¬ 
age  number  of  nodes  in  the  enumeration  trees  for 
both  algorithms,  when  using  relative  optimality  tol¬ 
erances  of  1%  and  5%.  For  Algorithm  Interdict,  this 
means  that  e  is  evaluated  as  0.01  [~ z(A(R),  R)"|  or 
0.05[~z(A(R), R)~\,  respectively,  instead  of  as  a  fixed 
value. 

These  tables  indicate  that  Algorithm  Interdict  runs 
between  one  and  two  orders  of  magnitude  faster  than 
the  hybrid  algorithm  (on  average,  24,  42,  and  119 
times  faster  on  the  grid  networks  Al,  A2,  and  A3, 
respectively).  We  note  that  Algorithm  Interdict  solves 
faster  than  reported  for  certain  problems  when  using 
other  schemes  for  determining  A  in  Step  1(f)  of  the 
Algorithm.  However,  we  have  not  optimized  algorith¬ 
mic  performance  by  specializing  the  code  for  particu¬ 
lar  problems.  Figure  2  illustrates  the  trade-off  between 


Table  2  Total  Time  to  Solve  BMXFI,  and  Average  Number  of  Nodes  in 
Enumeration  Trees  for  Grid  Networks  Al 


Value 

Tol.  (%) 

10x20 

Network  Al 

20x  40  30  x  60 

40x80 

Hybrid  alg. 

Tot.  time  (sec.) 

1 

1.2 

4.4 

59.7 

317.8 

Tot.  time  (sec.) 

5 

1.2 

4.4 

59.7 

319.8 

Avg.  nodes 

1 

0.5 

0.0 

0.1 

0.0 

Avg.  nodes 

5 

0.5 

0.0 

0.1 

0.0 

Alg.  interdict 

Tot.  time  (sec.) 

1 

0.1 

0.3 

2.6 

6.4 

Tot.  time  (sec.) 

5 

0.1 

0.7 

2.7 

6.5 

Avg.  nodes 

1 

6.1 

0.2 

0.4 

7.4 

Avg.  nodes 

5 

3.1 

0.2 

0.4 

3.4 

Table  3  Total  Time  to  Solve  BMXFI,  and  Average  Number  of  Nodes  in 
Enumeration  Trees  for  Grid  Networks  A2 


Network  A2 


Value 

Tol.  (%] 

1  10x20 

20x40 

30x60 

40x80 

Hybrid  alg. 

Tot.  time  (sec.) 

1 

1.2 

3.4 

459.8 

1,619.6 

Tot.  time  (sec.) 

5 

1.1 

3.4 

428.8 

1,053.0 

Avg.  nodes 

1 

0.1 

0.3 

0.5 

1.1 

Avg.  nodes 

5 

0.1 

0.3 

0.5 

1.1 

Alg.  interdict 

Tot.  time  (sec.) 

1 

0.1 

0.7 

5.9 

23.3 

Tot.  time  (sec.) 

5 

0.1 

0.7 

5.6 

13.3 

Avg.  nodes 

1 

13.5 

26.6 

56.6 

200.6 

Avg.  nodes 

5 

10.2 

5.6 

20.2 

31.7 

maximum  flow  and  total  interdiction  cost  for  the  case 
A3  (20  x  40).  The  lack  of  convexity  in  this  figure 
implies  that  the  efficient  frontier  cannot  be  identified 
though  weighted-sums  scalarization  alone  (Climaco 
et  al.  1997). 

The  second  collection  of  test  problems  uses  net¬ 
works  derived  from  the  highways  and  roads  in 
Maryland,  Virginia,  and  Washington,  D.C.  (Carlyle 
and  Wood  2005).  These  data  include  all  interstate, 
state,  and  county  roads  in  these  areas.  Using  vari¬ 
ous  arc  capacities  and  interdiction  costs,  we  construct 
14  different  network  problems  denoted  B1-B7  and 
C1-C7.  Each  problem  contains  multiple  sources  and 
sinks  that  are  connected  to  a  super-source  and  super¬ 
sink,  respectively,  using  invulnerable,  infinite-capacity 
arcs. 

The  "road-B  networks"  incorporate  9,876  directed 
road  segments,  which  are  modeled  as  arcs  in  a  di¬ 
rected  graph  with  3,670  nodes.  These  segments  rep¬ 
resent  all  roads  with  speed  limits  of  30  miles  per 
hour  or  higher.  16  nodes  along  the  region's  north¬ 
ern  border  constitute  source  nodes,  and  14  nodes  on 
the  southern  border  constitute  sink  nodes.  Includ¬ 
ing  auxiliary  nodes  and  arcs,  the  network  comprises 
3,672  nodes  and  9,906  arcs.  Interdiction  costs  and  arc 
capacities  are  generated  as  specified  by  Table  5,  where 
Prob(;i/;  =  /jl)  is  the  probability  that  uk  takes  on  a  value 
in  {1,  2, . . . ,  49).  Randomly  generated  data  are  statis¬ 
tically  independent.  Deterministic  problem  instances 


Table  4  Total  Time  to  Solve  BMXFI,  and  Average  Number  of  Nodes  in 
Enumeration  Trees  for  Grid  Networks  A3 


Value 

Tol.  (%) 

10x20 

Network  A3 

20x  40  30x  60 

40x80 

Hybrid  alg. 

Tot.  time  (sec.) 

1 

2.5 

196.5 

1,556.5 

2,257.0 

Tot.  time  (sec.) 

5 

2.5 

177.7 

651.1 

1,219.7 

Avg.  nodes 

1 

0.6 

1.4 

2.0 

0.4 

Avg.  nodes 

5 

0.6 

1.0 

0.7 

0.6 

Alg.  interdict 

Tot.  time  (sec.) 

1 

0.1 

0.9 

7.1 

28.3 

Tot.  time  (sec.) 

5 

0.1 

0.9 

5.4 

18.7 

Avg.  nodes 

1 

26.8 

100.1 

86.6 

92.1 

Avg.  nodes 

5 

19.3 

60.0 

25.7 

5.9 
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Figure  2  Solving  BMXFI  for  Case  A3  (20  x  40):  The  Trade-Off  Between 
Maximum  Flow  and  Total  Interdiction  Cost 

("Deter.")  have  uk  =  sk/5  for  all  k,  where  sk  is  the 
speed  limit  in  miles  per  hour  for  road  segment  k,  sk  e 
{30, 45, 50, 55, 65}.  Also,  for  B 7  and  C 7,  "Deter."  indi¬ 
cates  that  rk  =  1  when  sk  =  30,  and  otherwise  rk  =  2. 

The  "road-C  networks"  have  the  same  topology  as 
the  B-networks,  but  with  different  sources  and  sinks. 
152  nodes  in  the  region's  center  constitute  source 
nodes,  and  35  nodes  scattered  about  the  periphery 
constitute  sink  nodes.  The  network  comprises  3,672 
nodes  and  10,031  arcs  in  total.  Interdiction  costs  are 
generated  as  specified  in  Table  5. 

A  node  i  e  N  in  the  road  networks  is  said  to  be 
a  nonintersection  if  FS(i)  =  {(i,  ji),  (i,;2)}  and  RS(i)  = 
!(/i,0,  (/2/0l-  About  40%  of  the  nodes  in  the  road 
networks  are  nonintersections  because  of  the  level  of 
detail  represented:  Nonintersections  do  actually  rep¬ 
resent  intersections,  but  at  least  one  intersecting  road 
segment  has  been  deleted  because  its  speed  limit 
lies  below  the  cutoff  of  30  miles  per  hour.  Clearly,  a 
nonintersection  can  be  deleted  and  replaced  by  two 
arcs,  (/, ,  j2)  and  (j2 ,  /, ),  with  appropriate  adjustments 
to  the  numerical  data.  Both  algorithms  use  the  pre- 
processed  data.  (Solution  times  reported  exclude  the 
pre-processing  time,  but  none  exceeds  0.1  seconds.) 

For  the  cases  in  which  arc  capacities  uk  are  defined 
as  one  fifth  of  the  speed  limit  (see  Table  5),  Algo¬ 
rithm  Interdict  adds  a  test  to  ensure  that  the  lower 


Table  5  Parameters  for  Generating  Road-Network  Capacities  and 
Interdiction  Costs 


Network 

B1 ,  Cl 

B2,  C2 

B3,  C3 

B4,  C4 

B5,  C5 

B6,  C6 

B7,  C7 

Prob(r*  =  1) 

1 

1/2 

1/3 

1 

1/2 

1/3 

Prob(r*  =  2) 

0 

1/2 

1/3 

0 

1/2 

1/3 

Deter. 

Prob(r*  =  3) 

0 

0 

1/3 

0 

0 

1/3 

Prob(i4  =  ji) 

1/49 

1/49 

1/49 

Deter. 

Deter. 

Deter. 

Deter. 

Table  6  Total  Time  to  Solve  BMXFI,  and  Number  of  Nodes  in 
Enumeration  Trees  for  Road-B  Networks 


Value 

Tol.  (% 

)  B1 

B2 

B3 

Network 

B4  B5 

B6 

B7 

Hybrid  alg. 

Tot.  time  (sec.) 

1 

6.9 

12.0 

8.4 

6.9 

23.2 

13.3 

14.6 

Tot.  time  (sec.) 

5 

6.9 

12.0 

8.3 

6.8 

23.1 

13.1 

14.4 

Avg.  nodes 

1 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

Avg.  nodes 

5 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

Alg.  interdict 

Tot.  time  (sec.) 

1 

0.2 

0.3 

0.3 

0.4 

0.5 

0.6 

1.2 

Tot.  time  (sec.) 

5 

0.2 

0.3 

0.3 

0.4 

0.4 

0.6 

0.8 

Avg.  nodes 

1 

0.0 

0.1 

0.2 

0.1 

0.3 

24.1 

55.6 

Avg.  nodes 

5 

0.0 

0.1 

0.2 

0.1 

0.3 

24.1 

0.5 

Table  7 

Total  Time  to  Solve  BMXFI,  and  Average  Number  of  Nodes  in 
Enumeration  Trees  for  Road-C  Networks 

Value 

Tol.  (% 

Network 

)  Cl 

C2 

C3 

C4 

C5 

C6 

C7 

Hybrid  alg. 

Tot.  time  (sec.) 

1 

9.5 

13.4 

14.2 

7.9 

9.9 

14.9 

19.8 

Tot.  time  (sec.) 

5 

9.4 

13.2 

14.1 

7.9 

9.8 

14.7 

18.2 

Avg.  nodes 

1 

0.0 

0.2 

0.0 

0.0 

0.0 

0.0 

0.0 

Avg.  nodes 

5 

0.0 

0.2 

0.0 

0.0 

0.0 

0.0 

0.0 

Alg.  interdict 

Tot.  time  (sec.) 

1 

0.5 

0.7 

8.5 

0.5 

0.5 

2.9 

5.4 

Tot.  time  (sec.) 

5 

0.5 

0.7 

0.7 

0.5 

0.5 

2.4 

0.9 

Avg.  nodes 

1 

0.1 

33.8 

958.6 

0.1 

0.5 

278.7 

389.5 

Avg.  nodes 

5 

0.1 

33.8 

2.0 

0.1 

0.2 

220.5 

0.3 

bound  does  not  take  on  a  value  of  flow  that  cannot  be 
maximum,  i.e.,  z(A (R),  R)  £  {1, 2, 3, 4, 5,  7, 8, 14}.  For 
instance,  a  nominal  lower  bound  of  4  increases  to  6. 

Tables  6  and  7  display  solution  times  for  both  algo¬ 
rithms  applied  to  the  road-B  and  road-C  networks, 
respectively.  On  average.  Algorithm  Interdict  remains 
an  order  of  magnitude  faster  than  the  hybrid  algo¬ 
rithm  (30  and  15  times  faster  on  the  B  and  C  networks, 
respectively).  Flowever,  the  road-C  problems  include 
three  cases,  C3,  C6,  and  C7,  with  relatively  small 
speed-ups  at  the  1%-tolerance  level.  The  hybrid  algo¬ 
rithm  appears  only  modestly  slower  in  these  cases 
because  Algorithm  Interdict  requires  substantial  enu¬ 
meration,  while  CPLEX  requires  almost  none. 

On  average.  Algorithm  Interdict  requires  more  enu¬ 
meration  than  does  the  hybrid  algorithm,  no  doubt 
because  the  Lagrangian  lower  bounds  are  weaker 
than  CPLEX's  cut-enhanced,  LP-based  bounds.  An 
increased  tolerance  can  significantly  reduce  the  num¬ 
ber  of  branch-and-bound  nodes  for  Algorithm  Inter¬ 
dict,  while  the  number  of  nodes  for  CPLEX  stays 
small,  almost  independent  of  the  optimality  tolerance. 

7.  Conclusions 

This  paper  describes  a  new  procedure,  "Algorithm 
Interdict,"  for  solving  the  bi-objective  maximum-flow 
interdiction  problem.  The  algorithm  first  identifies  a 
large  portion  of  the  efficient  frontier  using  weighted- 
sums  scalarization  of  the  two  objectives  to  be  min¬ 
imized,  maximum  flow  and  total  interdiction  cost. 
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We  interpret  this  through  the  theory  of  Lagrangian 
relaxation.  A  specialized  branch-and-bound  proce¬ 
dure,  involving  partial  cut  enumeration,  then  identi¬ 
fies  any  missing  parts  of  that  frontier. 

We  have  compared  Algorithm  Interdict  to  a  hybrid 
algorithm  in  computational  tests;  the  hybrid  algo¬ 
rithm  calls  a  standard  integer-programming  solver 
instead  of  the  specialized  branch-and-bound  proce¬ 
dure,  when  the  Lagrangian  solution  from  Algorithm 
Interdict  is  not  e-optimal.  With  rare  exceptions,  our 
algorithm  is  one  to  two  orders  of  magnitude  faster 
than  the  hybrid  algorithm;  on  average  it  is  40  times 
faster.  Its  efficiency  results  from  the  fact  that  bounds 
and  s-t  cuts  are  computed  via  the  solutions  of  inter¬ 
related,  and  easily  solved  maximum-flow  problems. 
The  algorithm  may  require  more  enumeration  than 
does  linear-programming-based  branch  and  bound, 
but  that  enumeration  is  highly  efficient.  Algorithm 
Interdict  also  provides  the  benefit  of  not  requiring  a 
licensed  solver. 

We  might  attempt  to  reduce  enumeration  by  im¬ 
proving  the  Lagrangian  lower  bound.  Indeed,  it  is 
clear  that  the  standard  solver  achieves  better  bounds 
through  the  use  of  integer  cutting  planes.  If  such  cut¬ 
ting  planes  could  be  identified  and  Lagrangianized 
with  appropriate  multipliers,  this  could  improve  the 
lower  bound.  Wood  (1993)  identifies  some  problem- 
specific  cutting  planes  that  could  be  explored  for  this 
purpose. 
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